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ABSTRACT 


The aerobang maneuver, one of three types of aero-assisted 
orbital change methods, holds the possibility of reducing fuel 
consumption for orbital craft capable of atmospheric entry. 
It has been previously shown that different types of vehicles 
provide varying results over a constant heating rate 
trajectory. Further investigation into the optimization of 
the aerobang maneuver in this thesis includes the effects of 
using thrust vector control, the examination of the effects of 
increasing fuel mass fraction to increase orbital inclination 
changes, and the effects of that increase on both angle of 
attack and heating rate. The aerobang maneuver is shown to be 
capable of significant changes in orbital inclination in 
either a fixed heating rate or a fixed angle of attack mode 


for the Maneuverable Reentry Research Vehicle. 


ili 


Tore. 


IV. 


TABLE OF CONTENTS 


INTRODUCTION 


A. BACKGROUND 


B. TYPES OF SYNERGETIC ORBITAL MANEUVERS 


l. Aeroglide 
2. Aerocruise 
3. Aerobang 
C. SCOPESOESTHESTS 
D. FLIGHT VEHICLE 


E. APPLICATIONS 


EQUATIONS OF MOTION 
A. SPACECRAFT POSITION 


B. ORBITAL INCLINATION 


HEATING RATE AND CONTROL LAWS 

A. ATMOSPHERIC MODEL 

B. HEATING RATE MODEL 

C. DETERMINING THE ANGLE OF ATTACK 


D. THE HEATING RATE PROFILE 


SIMULATION PROGRAM . . e s s e s où 


A. THE MAIN PROGRAM 


iv 


19 
17 


14 


16 
16 
17 
18 
19 


21 


21 


B. THE SUBROUTINES 


lo 6 
3. The 
4. The 


5. The 


CNTRL Subroutine 

ACEL Subroutine 

ORB Subroutine 

WRT and HDR Subroutines 


RK4 Subroutine 


C. PROGRAM VERIFICATION 


+ 


We AEROBANG PERFORMANCE ANALYSIS 


A. THRUST 


VECTOR CONTROL 


1. Subcircular Case 


2. Supercircular Case 


B. INCREASING FUEL MASS FRACTION 


1. Subcircular Case 


a. 


B. 


Constant Heating Rate 


Constant Angle Of Attack 


2. Supercircular Case 


a. 


b: 


Constant Heating Rate 


Constant Angle of Attack 


VI. CONCLUSIONS AND RECOMMENDATIONS 


. . . . هم ه هو Ae CON@BUSTONS‏ 


« « ه €« هو و . ۰ B'SRECOMMENDATIONS‏ 


APPENDIX A PROGRAM ORBIT. . . 2 22. 


SAMPLE AERO.DAT INPUT FILE. . . 


22 
22 
22 
25 
25 
29 


26 


2 
28 
28 
30 
S 
Si 
31 
35 
40 
40 


42 


47 
47 


48 


‚50 


59 


SAMPLE DATA.DAT INPUT FILE. 


SAMPLE OUT.DAT OUTPUT FILE. 


LIST OF REFERENCES 


INITIAL DISTRIBUTION LIST 


vi 


.60 


-61 


62 


63 


Figure 


Figure 


Figure 


Figure 


Figure 


Figure 


Figure 


Figure 


Figure 


Figure 


Figure 


10 


11 


LIST OF FIGURES 


Aerodynamically Assisted Orbital Transfer 


The MRRV 


C,/C, and Angle of Attack for the MRRV 


The Coordinate System 


Spacecraft Force Diagram 


Aerobang and Aerocruise Efficiency 


Velocities and Heating Rates 


Flow Chart for the Main Program 


Flow Chart for Subroutine CNTRL 


Thrust Vector Angle Effects-Subcircular Case 


Thrust Vector Angle Effects-Supercircular Case 


Inclination For Constant Heating Rate 


XII 


10 


13 


20 


23 


24 


30 


33 


Figure 


Figure 


Figure 


Figure 


Figure 


Figure 


Figure 


Figure 


Figure 


Figure 


Figure 


Figure 


12 


13 


14 


15 


16 


17 


18 


19 


20 


21 


22 


23 


Magnification of the Constant Heating Rate 


Inclination 


Angles of Attack For Constant Heating Rate 


Inclination For 30° Angle of Attack 


Inclination For 35° Angle of Attack 


Inclination For 40° Angle of Attack 


Heating Rate For 30° Angle of Attack 


Heating Rate For 35° Angle of Attack 


Heating Rate For 40° Angle of Attack 


Inclination for Constant Heating Rate 


Angles of Attack For Constant Heating Rate 


Inclination for 30° Angle of Attack 


Inclination For 35° Angle of Attack 


NE d: 


33 


34 


36 


37 


37 


38 


39 


39 


41 


41 


42 


43 


Figure 


Figure 


Figure 


Figure 


24 Inclination For 40° Angle of Attack 


25 Heating Rates For 30° Angle of Attack 


26 Heating Rates For 35° Angle of Attack 


27 Heating Rates For 40° Angle of Attack 


ix 


43 


45 


45 


46 


LIST OF VARIABLES 


Acceleration of the flight vehicle. 


Vehicle acceleration in the normal (radial) 
direction. 


Vehicle acceleration in the tangential direction. 
Vehicle acceleration in the binormal direction. 
Coefficient of Drag. 

Coefficient of Lift. 


Zeroth order coefficient of binomial expansion for 
C 


ሺ 
First order coefficient for С 
Second order coefficient for E 


Zeroth order coefficient of binomial expansion for 
C 


iE 
First 0۲16۳ ۶ ori m 

Second order coefficient for C,. 

Drag on flight vehicle. 

Sum of all forces on the flight vehicle. 


Forces on vehicle along x’ axis in vehicle frame of 
ር. ህጅ 


Forces on vehicle along y’ axis in vehicle frame of 
rere rences 


Forces on vehicle along z’ axis in vehicle frame of 
reference. 


Local ۳ ጻነ 1 


Angular momentum of orbiting vehicle. 


X,Y,Z 


Y 2° 


Po 


Orbit ımelinatıon 


Unit vector along 2 axis of the inertial coordinate 
systen. 


Coefficient for heating rate formula. 
Lift on flight vehicle. 
Mass of the flight vehicle. 


Power coefficient for velocity in heating rate 
formula. 


Power coefficient for density in heating rate 
formula. 


Stagnation point heating rate of flight vehicle 
Geocentric radius of flight vehicle. 

Reference geocentric radius for atmospheric model. 
Reference surface area of vehicle. 

Eust of flight vehicle. 

Velocitv of flight vehicle. 

Cartesian coordinates for inertial reference frame. 
Cartesian coordinates for rotating reference frame. 
Cartesian coordinates for vehicle reference frame. 
Angle of Attack of flight vehicle. 


Coefficient for scale height in atmospheric 
model. 


Flight path angle, referenced to local horizon. 
Thrust vector angle. 

Right ascension of spacecraft. 

Densitv of the atmosphere. 


Reference densitv of the atmosphere. 


xi 


Ps 


Reference density of the atmosphere for atmospheric 
model. 


Angle of bank of the flight vehicle. 
Declination of spacecraft. 


Heading of flight vehicle based on angle from 
equatorial plane. 


Angular velocity of the vehicle's orbit. 


XI 


I. INTRODUCTION 


A. BACKGROUND 


E allorpiting spacecraft have the ability to correct 
and change their orbits, either within a given orbital plane 
or by transferring to other planes for better positioning 
Anding On the purpose of the spacecraft. These maneuvers 
are performed using the spacecraft’s propulsion and attitude 
control systems and are limited by the amount of fuel the 
Spacecraft can carry onboard. 

Executing a non-coplanar transfer requires more fuel than 
an orbital change within the original plane, such as 
increasing or decreasing orbital velocity to change orbital 
radius or altitude. In addition to velocity changes, a non- 
coplanar transfer requires a change in the direction of travel 
of the spacecraft. Missions requiring numerous plane changes 
are greatly limited by the fuel capacity of the spacecraft. 

Since the early 1960's, the possibility of using 
aerodynamic forces as a means of controlling both velocity and 
direction of a spacecraft has been studied [Ref. 1]. With a 
properly designed spacecraft, the Earth's atmosphere (or that 
of any other planet) could be used to assist the propulsive 


force of the spacecraft. Such a spacecraft would need to be 


able to use the lift forces generated to overcome the drag 
present to maintain flight and be constructed to withstand the 
surface heating effects that would occur within the atmosphere 
at near orbital speeds. Figure 1 shows a simple schematic 
drawing of what an aerodynamically assisted orbital transfer 


would entail. 


EXO-ATM. 
ORBIT.” 





Figure 1 Aerodynamically Assisted Orbital Transfer 


B. TYPES OF SYNERGETIC ORBITAL MANEUVERS 
Synergetic maneuvers are defined as maneuvers that use 


both atmospheric (aerodynamic) and propulsive forces. For 


orbital transfers between orbital planes, they have been 
divided into three categories: aeroglide, aerocruise, and 
aerobang. 
1. Aeroglide 

As the name implies, the aeroglide maneuver utilizes 
the aerodynamic forces present to create a gliding, unpowered 
prugectory. While fuel is required to deorbit prior to and 
reorbit after the maneuver, the aeroglide itself relies only 
Gamene Lift generated by the spacecraft's interaction with the 
atmosphere. As a consequence of this, the glide must not only 
be performed at a sufficiently low altitude to take advantage 
of the more dense atmosphere, but also at velocities high 
enough to maintain flight. This results in extremely high 
surface heating rates, perhaps beyond the capacity of the 
spacecraft's structure to absorb or dissipate [Ref. 2]. This 
disadvantage makes the aeroglide maneuver less appealing for 
nen coplanar orbital transfers by itself, although it could 
conceivably be used in conjunction with another technique. 

2. Aerocruise 

The aerocruise maneuver uses the spacecraft’s 
propulsive force as well as the aerodynamic forces present in 
a glide. The drag encountered by the spacecraft is exactly 
compensated for by the thrust generated by the spacecraft’s 
engines. This allows the spacecraft to maintain a constant 


heating rate and avoid the overheating problems present with 


the aeroglide. This requires both constant speed and 
altitude. The engine is throttled to maintain thrust 
sufficient to counteract drag, and the bank angle is adjusted 
to have a portion of the lift oppose the centrifugal force, 
thus maintaining a constant orbital radius or altitude. 
Changes in orbital inclination come as a result of the lift 
generated and will depend on the lift to drag ratio (L/D) and 
the angle of attack (a) of the vehicle. Altering these 
parameters will produce varying amounts of inclination change 
for a given amount of fuel burned. 
3. Aerobang 

This maneuver also uses a constant heating rate as a 
requirement, allowing the altitude and velocity to vary eum 
the flight as dictated by the control law for the maneuver, 
discussed later. This allows a constant thrust to be used and 
allows the angle of attack to vary. The thrust is set at the 
maximum, reducing the time the spacecraft is in the atmosphere 
compared to the aerocruise, thus reducing the integrated 
heating effects on the spacecraft. It has previously been 
shown that, for certain velocities, altitudes, and heating 
rates, an aerobang maneuver will yield a greater orbital 
inclination change than an aerocruise maneuver, for the same 


fuel expenditure [Ref. 3]. 


C. SCOPE OF THESIS 

This thesis will further examine the aerobang non-coplanar 
orbital transfer maneuver. The initial results showed that 
the maneuver generated smaller than expected inclination 
changes for the given two percent mass fraction of fuel 
consumed. The thrust acted in the direction of the 
spacecraft’s body and thus only a portion (depending upon the 
angle of attack) acted to help the lift, which produced the 
inclination change. The effects of thrust vector control 
(TVC) are included here to allow the thrust to contribute 
further to the lift, and thus increase the inclination change 
while maintaining an angle of attack within the vehicle's 
capability. 

The aerobang maneuver is limited by the amount of fuel 
allotted for each atmospheric pass. Obviously, the greater 
the fuel used, the longer the maneuver can be performed, 
subject to heating and altitude restraints. The initial 
results [Ref. 3] were based upon a 0.02 mass fraction because 
the aerobang maneuver was unable to maintain a constant 
heating rate profile beyond that point for some conditions of 
altitude and velocity. The inclination changes generated with 
this restriction (less than one degree for the most part) is 
not useful for many orbiting spacecraft applications. By 
expanding the allowable mass fraction of the fuel consumed for 
the maneuver to values up to twenty percent, larger, more 


useful inclination changes should be able to be produced. 


In addition to the two limiting heating rates examined in 
Reference 3, the effects of both aerocruise and aerobang over 
the entire range of heating rates are studied. The intent is 
to determine the shape of the heating rate profile for both 
maneuvers and determine more precisely when aerobang is a 
superior method of orbital transfer. This study will also 
observe the behavior of the vehicle’s angle of attack during 
extended constant heating rate maneuvers as well as fixing the 
angle of attack for other maneuvers and observing its effect 


on the heating rate. 


D. FLIGHT VEHICLE 

The vehicle used for this study is the Maneuverable 
Reentry Research Vehicle (MRRV). The MRRV was designed as a 
hypersonic research vehicle to be launched by the space 
shuttle and capable of a powered re-entry and flight through 
the atmosphere to demonstrate the feasibility of aerocruise 
maneuvers (Ref. 4]. The MRRV has also been used as a test 
vehicle for aerobang simulations [Ref. 3]. Because of this, 
and for consistency, it is used here as well. 

As depicted in Figure 2 [Ref. 3], the MRRV is a winged 
body approximately 7.6 meters long with a four meter wingspan. 
The effective surface area is 11.61 m° with an initial mass of 


4898 kg, of which 2588 kg is fuel. The MRRV is powered by 


three Marquardt R-40-B rocket motors that provide a total 
thrust of 14679 Newtons with a specific impulse (1.0) of 295 


seconds. 





Figure 2 The MRRV 


Figure 3 shows that the MRRV has a nominal lift to drag 
ratio (L/D) of 2.3 at its optimal angle of attack of 14 
degrees. This is based upon wind tunnel data gathered out to 
an angle of attack of between 30 and 40 degrees [Ref. 4]. No 
data was available for higher values, implying that this is 


the upper limit for angle of attack for this vehicle. 


10 20 30 


Angle of Attack (Deg) 





Figure 3 C,/C, and Angle of Attack for 
the MRRV 


E. APPLICATIONS 

Initially, it would appear that the results described in 
Reference 3 (and modified later in this thesis for the 2% 
case) may not be useful for most orbital transfer applications 
Since they amount to less than a one degree change in 
inclination, occurring over a period of under 20 seconds. 
There are, however, possible applications and factors to take 
into account when studying these results. 

An inclination change of even a degree or so can 
correspond to a difference of over a hundred kilometers when 
the angle is projected up to the orbital altitudes used here 
and higher as well. For example, for a one degree change in 
inclination at a radius of 6447 km, the spacecraft will move 
112.5 km. In a military application involving pursued and 
pursuing vehicles in orbit, such as anti-satellite warfare, a 


quick aerobang maneuver resulting in even a small amount of 


orbital transfer may be capable of preventing a vehicle from 
being detected or destroyed. 

Even for civilian applications, a small aerobang maneuver 
may be preferred to a longer one to minimize the heat load 
imparted to a vehicle that may not have the thermal protection 
system to adequately protect against a longer atmospheric 
flight that would yield a greater inclination change. By 
making multiple, low fuel passes and reorbiting between 
orbital changes to radiate the heat to space, the same total 


inclination change can be made if time is not of the essence. 


II. EQUATIONS OF MOTION 


The motion of the aerodynamic spacecraft within its orbit 
is defined by six variables for both position and velocity. 
Position is defined by the standard spherical coordinates; r 
(radius), 8 (right ascension), and ¢ (declination). The other 
three variables are v (velocity), y (flight path angle) NEMO 
V (heading angle). These variables are graphically shown in 


Figure 4 below. 





Figure 4 The Coordinate System 
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ሮስ 1 path angle is described as the angle of the 
velocity vector with relation to the local horizontal plane 
tangent to the orbit and the heading angle is the direction 
angle within that plane relative to the local inertial 
latitude.These equations were derived (Ref. 3) assuming a 
spherical, non-rotating Earth. Atmospheric motion (winds) and 
the oblateness of the Earth were also neglected. The method 
of derivation and the resulting equations that were used in 


this studv are summarized in this chapter. 


A. SPACECRAFT POSITION 

The derivation of the equations of motion requires the use 
of three coordinate systens. As shown in Figure 4, the 
position of the vehicle can be expressed in the inertial XYZ 
Men pene orbital xyz system, also centered at the Earth's 
center, with the positive x axis pointed to the orbiting 
spacecraft, as well as the x’y’z’ system which is centered on 
the spacecraft with the positive y’ direction in line with the 
metocitv vector and the x' axis within the xy plane. 

Direction cosine matrices (DCM) are required to transform 
the various coordinate systems into one to be used for the 
expression of the spacecraft's position and velocitv. In 
order to get from the inertial XVZ svstem to the orbital xvz 
system, a bodv 3-2 rotation through angles 0 and ġ is 
required. iG @btain the xyz coordinates from the x'y'z' 


system a 3-1 rotation is used through angles y and y. 


11 


By equating the velocity vector to the time derivative of 


the position vector for the spacecraft, the kinematics can be 


expressed as follows: 


dr 


Zn (25 
46. Vcosvcosw (2.2) 
dt rcosġ 
dò _ Vcosysiny 209 
dt T 


The other three equations are based on the balance of 
forces on the spacecraft, and their relation to the 
accelerations and masses (F=ma). The forces are first derived 
in the spacecraft’s frame of reference, the x’y’z’ system, as 


shown in Figure 5. 


The forces in the x', y’ and z’ directions are: 


F,=[Tsin(a+e) +L)coso (2.4) 
ድ መሮን ከ و‎ (2.5) 
۳,2۶] 152112 )0 ۲6۱1 (226) 


Again, the conversion to the rotating reference, the xyz 
system, is done using the appropriate DCM. The resulting 


expressions for the forces on the spacecraft are then set 


12 





Figure 5 Spacecraft Force Diagram 


equal to the product of the mass and the accelerations 
(calculated by taking the derivative of velocity vector) in 
each of the three direction components (x, y, and z) and are 


meduced to the following: 
[Tcos (a *e) -D] -m,gsiny-m av (287) 


zo: 


dy_ [Tsin(a+e)+L]coso-m,gcosy, v? 


አፈ 2.8 

= m 7 COSY ( ) 
۰ Q 2 

y. di . [Tsin(a*e) *L] sino | V eo e cosytand (2.9) 
dt m,cosy I 


The components of the spacecraft's acceleration in the 
tangential, normal, and binormal directions can be expressed 


independently. This will allow the final three equations of 


13 


motion to be expressed in a more compact form. For Ма 
perfectly, circular orbit, the normal (to “the tang. MEE 
velocity and acceleration vectors) direction is also the 
radial direction, pointing to the center of the Earth. Hence, 


the notation A, for this component. 


_ [Tsin(a+e) +L]coso 


A ۰ 
R m, (2.10) 
A ee (2258 
Ma 
A O (2.12) 
Ma 
The final three equations of motion then become: 
dv ; 
—— =A.-gsin 2 208 
ን 95 7۷ ( ) 
d Ve 
ር ae (2.14) 
y AY E S ict (25256 
dti 9۹5۳۲ E 


B. ORBITAL INCLINATION 

The measure of the effectiveness of the orbital plane 
change maneuver is the change in the orbit’s inclination. 
Since the equations of motion use the spherical coordinates to 
describe the spacecraft's motion, the inclination must also be 


expressed in these coordinates. As shown in Reference 3, the 


14 


ongulorimomentumivector will be perpendicular to the orbital 
plane so that the angle between that vector and the unit 
Ше сог KE ЖЕН ПОЛ the positive 4 direction will be the 
orbit's inclination. The magnitude of the angular momentum 


Mector Ys "shown to be Vrcosy while the Z component of the 


angular momentum is Vrcosycosy. Since 
b'kdHcosi (2.16) 
chen 
Vreosycoswcosd=|Alcosi (2509) 


and so it follows that 


cosi=cosycosb (2.18) 


giving a means of simply calculating the orbital inclination 


for any plane change maneuver. 
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III. HEATING RATE AND CONTROL LAWS 


The control parameters for the aerobang maneuver are the 
bank angle, the angle of attack, and the thrust. The angle of 
bank is set to 90 degrees, the thrust is set to the maximum, 
and the angle of attack is adjusted to maintain a constant 
heating rate on the spacecraft during the maneuver. In order 
to implement a control law for the aerobang orbital transfer 
maneuver, the modeling of the atmosphere and resulting heating 


rate must first be developed. 


A. ATMOSPHERIC MODEL 

The first step in the modeling of the heating Tate 
behavior is to model the atmosphere in the region in which the 
maneuvers will take place. For the aerobang maneuvers 
examined in this and previous studies, the region between the 
altitudes of 50 and 120 km was used. The density of the air 
at these altitudes was modeled using the U.S. Standard 
Atmosphere, 1976, [Ref. 5, p. 7] and is approximated by a 


local exponential model 


P 079 (3000‏ و م-م 


where 


p. = 3.0968x10” kg/m’ 
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በ96 een m! 


r. = 6438 km 


5 


B. HEATING RATE MODEL 

Atmospheric density (and therefore altitude, as shown 
above) along with the vehicle's velocity combine to give a 
basic relationship for the heating rate. As shown in 


Bererence 5, that relationship is: 

0-۳ (372) 
with k, n, and m as constants. This equation, combined with 
the equations of motion can be used to derive the angle of 
attack required to fly at a constant heating rate during the 
aerobang maneuver. 

By combining Equations 3.1 and 3.2, the relationship 
between density, heating rate, and velocity becomes more 
apparent. 


۵ 


. ር [m (ars) 
I 


Ve 








In order to show the effect of the heating rate model on the 
equations of motion, the derivative of Equation 3.3 is taken, 


and simplified to give: 
ALA nat 


Using Equation 2.1, the expression for dr/dt can be used to 


EXPRESSO V/at as a function of V and y, the flight path angle. 
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dV 2 | دوس‎ 
2 9፡7 sinv (3:59 


Using Equations 2.11 and 2.13 gives for the equation along the 


path, 


Tcos (a+e) -D 


Ma 


dv ۱ 
SaS -gsin 3 0 
gsiny (3.6) 
C. DETERMINING THE ANGLE OF ATTACK 

Equating Equations 3.5 and 3.6 gives an expression for 
calculating the angle of attack required for the aerobang 


maneuver: 
0=Tcos (a+e) -D-m,siny (g+ £@ v2) (399p 


The drag term, D, is also an function of a so the solution for 
alpha is given by the transcendental equation above. 

In order to solve for an angle of attack to BE sS 
throughout the maneuver, an initial value must be selected and 
the above equation (3.7) must be iterated to solve for a at 
each successive time interval. As shown in Reference 3, 
Newton's Method of convergence was used to arrive at an angle 
of attack for this situation. Initially guessing a involved 
using a binomial curve fit for the drag versus angle of attack 
data. Assuming minimal values for the tangential acceleration 


and angle of attack (so that the cosine is nearly equal to 
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er sendenoschrust vector angle in Equation 2.11, then the 
thrust is equal to the drag (or the binomial expansion curve 


INE for drag), 


Дд 


ш = 
15 DIAS 


(3.8) 


where C,, C,, and C, are the coefficients for the zeroth, first 


57 
and second order terms, respectively, in the binomial curve 
fit expression (Ref. 3: p. 13]. The quadratic formula gives 


a method for calculating the first guess for a. 


(559) 
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This technique is used in the computer program to simulate 
the aerobang maneuver to get a satisfactory approximation for 


the varying angle of attack. 


D. THE HEATING RATE PROFILE 

The previous thesis examined the aerocruise and aerobang 
at only two heating rate values. It was assumed that the 
resulting regions of superior performance for both subcircular 
and supercircular speeds held for all heating rates in between 
let. 3: p.63]. 

Prior to a closer examination of the behavior of flight 
parameters for the aerobang maneuver, this assumption was 


validated by determining the velocities at which aerocruise 
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lost its superiority to aerobang in terms of the inclination 
changes produced during the 2% mass fraction runs. These 
results, shown in Figure 6, do show that at very near circular 
speeds, the aerocruise maneuver did produce better results, 
but that the profile is not entirely symmetrical about the 
circular velocity: In general it was true that, at lower 
heating rates, the aerocruise method was better for a much 


smaller range of velocities. 
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Figure 6 Aerobang and Aerocruise Efficiency 
Velocities and Heating Rates 
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IV. SIMULATION PROGRAM 


In order to evaluate the potential of the aerobang orbital 
plane change maneuver, a FORTRAN program, previously developed 
(Ref. 3), was modified and used to simulate the action of the 
spacecraft while in the atmosphere. The modifications to the 
original program included the addition of e, the thrust vector 
angle, the addition of a constant angle of attack option, and 
the elimination of the aerocruise option and the associated 
calculations. 

The program, listed in Appendix A, was written in a 
modular format, using several subroutines for calculations, 
integration of the equations of motion, and output. The main 
program calls each one and controls their use in the program. 
For the heating rate versus velocity simulations used in 
Chapter III, the original program [Ref 3; App. A] was used to 


produce the aerocruise data shown in Figure 6. 


A. THE MAIN PROGRAM 

The main portion of the program links the subroutines 
together and controls the iteration of the solution. The two 
input files and one output file (also in Appendix A) are 
opened first. The use of the two input files (AERO.DAT and 
DATA.DAT) minimized the changes that had to be made to the 


program and subroutines themselves. The AERO.DAT file 
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contains the information about the spacecraft and the DATA.DAT 
file contains the information concerning the flight profile 
and the printing interval. The main program also calculates 
the mass changes throughout the simulation and stops the 
program when the final mass is reached or if the time extends 
beyond the maximum allowed. A flow chart for the main program 


appears in Figure 7 [Ref. 37: 


B. THE SUBROUTINES 
1. The CNTRL Subroutine 
The control subroutine is the major subroutine, 
calculating the atmospheric and heating rate models, as well 
as the angle of attack for the aerobang control itself. The 
angle of attack is computed using the Newton approximation 
routine as discussed in Chapter II to simulate the modulation 
of the angle of attack for the  aerobang maneu eE: 
Alternatively, the angle of attack can be fixed and allowed to 
change the heating rate accordingly. In either case, the lift 
and drag forces are calculated here for use later in the 
equations of motion. Figure 8 shows the flow chart for this 
Subroutine. 
2. The ACEL Subroutine 
This subroutine computes the three accelerations, 
defined by Equations 2.10, 2.11, and 2.12. These values are 
used in the next subroutine to evaluate the motion of the 


5۵ 16 6 0 E. 
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Program orbit starı 


Open input files 
Vehicle aerodvnamics 
Orbital parameters 
Open output file 


Integrauon of orbital vehicles 
flight path for 1/4 time step Has vehicle exhasted 
fuel or reached end 
of time interval? 


Is Runge-Kutta 
routine done 
with integration? 





Figure 7 Flow Chart for the Main Program 
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Begin Subroutine 


Define Constants 


Calculate Local 


Initial Guess for 


Angle of Attack 


Calculate Lift and 


Calculate Change to 
Angle of Attack 


Is the Change 
Larger 


Changed? 


RETURN 


Figure 8 Flow Chart for Subroutine CNTRL 
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3. The ORB Subroutine 
In this portion of the program, the six equations of 
ከህ ፣ክቢ፡ቢንቢን ቢን) ዕቕይ፡ ያ ሽበን. 2.3, 2.13, 2.14, and 2.15) are 
defined in terms of the six variables defining the motion of 
the spacecraft; given as a vector, X(6), along with their 
derivatives, XDOT(6). These variables are used in the 
subroutine to compute the motion of the spacecraft at each 
time iteration. The equations are integrated in a Runge-Kutta 
routine, described later. 
4. The WRT and HDR Subroutines 
These subroutines control the output of the program. 
The HDR (header) routine creates a header in the output file 
that gives the initial conditions of the program and the 
headers for the columns of data that the program computes. 
The WRT (write) subroutine is called from the main program and 
used on each time iteration to write the calculated data into 
the output file. It is also called at the beginning and end 
Srzehesiterations to put both the initial and final data rows 
into the file. 
5. The RK4 Subroutine 
This subroutine is the fourth order Runge Kutta 
routine used to integrate the equations of motion. It was 
written by Professor I. M. Ross and is a part of the public 
collection of subroutines at the Naval Postgraduate School. 


Here, the time iterations are counted and used to create the 
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multiple integration computations and send the results to the 


main program for writing into the output file. 


C. PROGRAM VERIFICATION 

Since this program is based on an earlier version, the 
validation of this program involved using it to verify results 
obtained in the previous study. By setting the thrust vectör 
angle equal to zero in the AERO.DAT file's last entry, the 
exact same results should be obtained for similar initial 
conditions used previously. 

The comparison of the aerobang and aerocruise maneuvers to 
verify the applicability of aerobang over the original heating 
rate profile provided an opportunity to verify the program's 
proper operation. By using both programs to calculate Fes 
for the aerobang cases and obtaining the same results, the new 
program could then be used on its own to incorporate the 
effects of thrust vector control and serve as a tool for 


further analysis of the aerobang maneuver. 
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у. AEROBANG PERFORMANCE ANALYSIS 


As shown in Chapter III and previous studies [Ref. 3], 
Mere exist certain conditions in which the aerobang orbital 
plane change is the most efficient as determined by the amount 
of orbital inclination change that can be achieved for a given 
amount of fuel expended. Having determined the conditions 
necessary for successful aerobang maneuvers, the next step is 
to determine the ways in which the aerobang maneuver can be 
improved and best utilized to achieve a given mission 
eriterion. The effects of adding a thrust vector control 
(TVC) system to improve upon previously obtained optimization 
results was examined first. Further simulations focused on 
the effects of increasing the mass fraction of fuel burned to 
allow for more useful inclination changes (on the order of ten 
or more degrees per maneuver) and the advantages and 
disadvantages of fixing the angle of attack rather than the 
heating rate for a specific maneuver. 

In order to remain consistent with the previous study 
[Ref. 3] and to be able to compare the previous results with 
those obtained and given here, a review of those results is 
necessary. 

The previous thesis looked at both subcircular and 
supercircular velocities for both aerocruise and aerobang 


maneuvers for two distinct heating rates. The effects of the 
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heating rate values over the entire range rather than just the 
two endpoints on the viability of both subcircular and 
supercircular results were discussed in Chapter III. These 
results were based on a mass fraction of fuel expended of only 
tworperceent. 

For the higher heating rate (102255108 W/m^), the best 
aerobang performance results, in terms of inclination changes, 
occurred at a velocity factor (k-V /V) of 0.98 Eo ۳9 
supercircular case and 1.02 for the subcircular case, with 


inclination changes of 0.554? and 0.789?, respectively. 


A.  THRUST VECTOR CONTROL 
1. Subcircular Case 

In order to determine the overall effects of thrust 
vector control on the aerobang maneuver, the computer 
simulations were run for a range of e from -80 to 80 degrees. 
while these angles may not presently be feasible from a 
structural or control point of view, they were includes 
completeness. The simulations were run at an initial radial 
distance of 6447 km, corresponding to an altitude of 77 km, 
and a heating rate of 1.42 MW/m°. As shown in Figure 9 below, 
the most efficient value for & fell between 10? and 16?, giving 
an inclination change of 0.8? over that range of vectoring 
angles. For negative values of thrust vector angles, the 
inclination changes were not as large as those at zero and 


positive angles since negative values of e€ will have a 
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component in the direction directly opposing the lift which is 


responsible for the change in inclination. 


THRUST VECTOR ANGLE/INCLINATION SUBCIRCULAR CASE 
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Figure 9 Thrust Vector Angle Effects - Subcircular Case 
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2. Supercircular Case 

For k = 0.98, at the same altitude and heating rate, 
similar results were obtained for the addition of thrust 
vectoring for this case. As shown in Figure 10, the best 
inclination change occurred with a thrust vector angle of 15 
to 20 degrees, resulting in an inclination change of 0.570°. 
At higher e, the angles of attack decreased as the spacecraft 
moved faster to higher less dense air with a subsequent loss 


OE LIFE, 
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Figure 10 Thrust Vector Angle Effects - 
Supercircular Case 
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The shape of the curve over the broad range of thrust vector 
angles is again indicative of the effects of extreme values of 


e on flight performance. 


B. INCREASING FUEL MASS FRACTION 

۳۰ ۱9 Be ۱9۳/9 to conclude that by increasing the 
amount of fuel consumed for an aerobang maneuver, a greater 
change in inclination would be achieved. By establishing a 
goal of achieving an inclination change of 10° or more, the 
effects of increasing the allowable mass fraction were 
investigated. These simulations were run for both subcircular 
and supercircular orbits (k=1.02 and 0.98), at an initial 
altitude of 77 km (R=6447 km), and restricted to the same 
heating rate as was used in the previous 2% fuel mass fraction 


runs. 


1. Subcircular Case 
a. Constant Heating Rate 
Initially, the increase of mass fraction was tried 
for the constant heating rate case usually prescribed for the 
aerobang maneuver. As the mass fraction was increased from 
the original 2% up to 20%, it became apparent that the 
technique of modulating the angle of attack to maintain a 
constant heating rate profile was difficult to implement. 
Figure 11 shows the expected effect of increasing mass 


fraction on aerobang maneuvers. The inclination change is 
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very impressive, between 25° and 28°, depending upon the thrust 
vector angle. In these simulations, five values of e (0, 10, 
20, 30, and 40%) were used to show the effects of TVC on these 
extended runs. The thrust vector angle labels for this figure 
refer to the right hand side of the curves; zero at the top, 
40 at the bottom. The smaller values of e, gave the greater 
inclination changes. While this seems to contradict the 2% 
mass fraction results, the data for time values less than 60 
seconds revealed that the behavior of the inclination followed 
the behavior shown earlier for the 19 second, 2% mass fraction 
case. The e = 10° curve starts out above the 0° curve, but the 
curves cross at 56 seconds and the 0° curve finishes slightly 
above the 10° curve. The two curves are nearly 
indistinguishable from each other, but the results for the 2% 
case (the slight advantage of the ten degree thrust angle with 
the three other cases even lower) are shown in Figure 12 and 
are consistent with the results for the longer times of Figure 


PA 
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Figure 11 Inclination For Constant Heating Rate 
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Figure 12 Magnification of the Constant Heating Rate 
Inclination 
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This superior performance in orbit transfer comes with a 
penalty. Figure 13 shows the change in angle of attack over 
the extended simulations. Values of a rise to between 65 and 
70 degrees, depending again on the thrust vectoring anma 
decrease to about 53 or 54°. The zero thrust angle case had 
the highest angles of attack throughout the flight. The 
addition of thrust vectoring caused the angle of attack to 
decrease. This reduction in a reduced the lift generated and 


the inclination changes as well. 
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Figure 13 Angles of Attack For Constant Heating Rate 
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The aerodynamic data provided for the MRRV, mentioned in 
Chapter I, only covered angles of attack of up to 40°. It can 
be assumed, therefore, that the vehicle is probably not 
capable of flying at a’s greater than 40°. If that is the 
case, then a vehicle of the same surface area, initial mass, 
Mier Durning rate, and thrust would not meet these high angle 
Of attack results, if they are valid, while maintaining the 
constant heating rate profile simulated in this study. 

b. Constant Angle Of Attack 

АС point, a departure from the "standard 
aerobang maneuver theory was attempted to try to avoid the 
angle of attack problems described above. Rather than 
allowing the angle of attack to vary and allowing the constant 
heating rate control the maneuver, the angle of attack was 
fixed and the heating rate was allowed to vary to evaluate the 
possible inclination changes while observing the heating rate 
to see if it would exceed the desired heating rate and, if so, 
by how much. 

An inclination change of approximately 10° was still a 
goal and several values of a were tested to see if this goal 
could be met. The graphs in Figures 14, 15, and 16 show the 
results of these runs. Angles of attack of 30°, 35°, and 40°, 
respectively, were selected to get the desired ten degree 
inclination change. Again, the same five values of thrust 


vectoring are included to show their effects. As a increases, 
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the inclination changes rise from values of 6.5° to 9.5° at a 
= 30° to values of 12° to 17° at a = 40°. In this case, the 
increasing thrust angle with the fixed angle of attack 
resulted in increasing inclination changes throughout the 
runs. With a fixed, its contribution to the vehicle’s lift is 
fixed and the changes in inclination for this case are made by 
the contribution of the thrust vectoring to спе 8 


generated. 
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Figure 14 Inclination For 30° Angle of Attack 


36 


40 deg 
JO deg 


መጫ 
a 
ወ 
o 
| 
Z 
Q 
< 
zZ 
5 
e 
= 








TIME (sec) 





Figure 15 Inclination For 35° Angle of Attack 
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Figure 16 Inclination For 40° Angle of Attack 
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The result of more concern is that of the heating rates 
that existed as a result of the fixed a flight regimes 
Figures 17, 18, and 19 show the heating rates obtained for 
these simulations. They vary from 1.8 x 10% W/m? on up to 2.2 
x 10% W/m, both significantly above, almost up to twice thé 
value used for all previous runs. A vehicle used to perform 
these fixed a aerobang maneuvers would, therefore, need to 
have a thermal protection system capable of withstanding these 


heating rates. 
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Figure 17 Heating Rate For 30° Angle of Attack 
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HEATING RATE (W/M2 X E+07) 
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Figure 18 Heating Rate For 35° Angle of Attack 
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Figure 19 Heating Rate For 40° Angle of Attack 
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2. Supercircular Case 
a. Constant Heating Rate 

Supercircular aerobang maneuver simulations using 
20% fuel mass fraction were attempted over the same range of 
angles of attack and thrust vector angles as for the 
subcircular case. Under otherwise identical conditions, the 
supercircular simulations were unable to generate any results 
that were even close to those of the subcircular case. The 
inclination changes were no greater than 0.865° for an e of 10° 
and the burn time was no more than 37 seconds, approximately 
4% mass fraction as opposed to the desired 20% set in the 
simulation progran. The combination of decreasing a and 
increasing altitude at the supercircular velocity caused the 
vehicle to lose its lifting force, resulting in the shortened 
flight. Compared to the subcircular case, the supercircular 
aerobang maneuvers proved to be ineffective in producing 
substantial inclination changes (on the order of 10°) for any 
kind of significant increases in fuel mass fraction, although 
they may still be viable for short, quick maneuvers. Figures 
20 and 21 show the inclination and angle of attack results for 


this case. 
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Figure 20 Inclination for Constant Heating Rate 
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Figure 21 Angles of Attack For Constant Heating Rate 
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b. Constant Angle of Attack 
This case proved to be more successful since the 
angle of attack (and therefore the lift generäted) never 
decreased to the values near zero as it did in the constant 
heating rate case. Each aerobang simulation was completed for 
the entire 20% burn and significant inclination changes were 


observed. Figures 22, 23, and 24 show these results. 
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Figure 22 Inclination for 30° Angle of Attack 
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Figure 23 Inclination For 35° Angle of Attack 
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Figure 24 Inclination For 40° Angle of Attack 
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As in the subcircular case, the inclinations here 
increased with higher angle of attack and with increased 
thrust vector angles. The maximum values, however, were not 
as large since the supercircular trajectories caused the 
altitude of the vehicle to increase and the resulting 
atmospheric density to decrease, thus reducing the lift 
Capability. 

The heating rates for supercircular flights, however, 
Showed an inherent advantage with this case. The higher 
altitudes had a greater effect on heating rates than the 
increased velocities did, resulting in decreasing values for 
heating rate. As shown in Figures 25, 26, and 27; the heating 
rates for all cases started with the same value of 1.45 x 10° 
W/m°. The final values depend on both angle of attack and 
thrust vector angle. The higher angle of attack cases 
resulted in less pronounced decreases in heating rate Vale 
from start to finish. The thrust vector angles in these cases 
actually contributed to higher heating rates as well, due to 
the slight increase in lift and velocity that they provided. 
In particular, as seen in Figure 27, a limit is reached Tn WE 
40? thrust vector case, where the heating rate reaches a 


minimum limit and then begins to increase again. 
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Figure 26 Heating Rates For 35° Angle of Attack 
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Figure 27 Heating Rates For 40° Angle of Attack 
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VI. CONCLUSIONS AND RECOMMENDATIONS 


A. CONCLUSIONS 

The aerobang orbital plane change maneuver has been shown 
to have both commercial and military applications for 
obtaining significant inclination changes for aerodynamically 
assisted orbital transfers. The improved efficiency of the 
aerobang technique over the aerocruise for subcircular and 
supercircular speeds at two discrete heating rate values given 
in Reference 3 was shown here to extend to the entire range of 
heating rates studied, not just the high and low heating rate 
cases. 

Thrust vector control can be effective in increasing 
inclination changes for both short and long duration aerobang 
maneuvers. This was demonstrated for both supercircular and 
subcircular velocities. The use of thrust vector control was 
most effective when applied at a constant angle of attack. 

The effects of increasing mass fraction of fuel consumed 
from two to twenty percent had a significant effect on the 
aerobang maneuver, beyond the obvious extension of the time 
used to perform the maneuver. For the subcircular maneuvers, 
ee sconsrant heating rate method of controlling the aerobang 
flight resulted in excellent inclination changes but with 


extremely high (70°) angles of attack. By controlling the 
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angles of attack, the heating rates increase by nearly a 
factor of two, placing new requirements on any spacecraft 
design. 

The supercircular aerobang maneuvers were also simulated 
for both constant heating rate and constant angle of attack 
cases. The constant heating rate case proved to be 
ineffective in generating even one degree of inclination 
change because of the loss of angle of attack and, therefore, 
lift. The constant angle of attack runs, however, proved to 
be almost as effective as for the subcircular cases. This 
case also provided two inherent advantages over the 
subcircular case, an increasing altitude which translates to 
a decreasing heating rate value and the likelihood of an 


easier reorbit maneuver. 


B. RECOMMENDATIONS 

There remain several areas in the study of aerodynamically 
assisted orbital plane change maneuvers that have yet to be 
explored. Among these are the effects of combining any two or 
all of the three techniques (aeroglide, aerocruise, and 
aerobang) to maximize the fuel efficiency of a given 
inclination change requirement. 

This thesis, as well as most other work in this area, 
deals only with the atmospheric passage portion of the orbit 
transfer. By extending the simulations to cover an entire 


orbit transfer, including both entering and exiting the 
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ae pere, a Detter measure Of the effectiveness of these 
maneuvers can be made. 

This is the first work in which the constant angle of 
attack method of controlling the aerobang maneuver was used. 
Purther studies should concentrate on whether or not it can 
indeed be a viable means of control, from both a flight 
control system and thermal protection point of view. 

їп Order to verify the applicability of aerodynamically 
assisted maneuvers, the study of thermal protection system 
requirements and techniques must be included. Superior 
atmospheric and heating models can also be used to refine the 
requirements for such systems. 

Finally, a means of minimizing the time required for a 
maneuver to meet a given inclination change requirement could 
be established. This would allow the mission goal to be met 


Whale minimizing the integrated heat load on the spacecraft. 
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APPENDIX A 


PROGRAM ORBIT 
IMPLICIT REAL*8(A-H,L-Z) 
DIMENSION X(6),XDOT(6),C(6) 


PROGRAMMER: THOMAS P. SPRIESTERBACH 
MODIFIED BY: RICHARD E. JOHNSON 


DATE 
SYSTEM 


I 
INDEX 
J 
KOUNT 
DUM 
OLDCO 


CHANGE 
NE 


BE 


15 FEB 93 
IRIS INDIGO/FORTRAN 


VARIABLE LIST 


COUNTING INDEX 

COUNTING INDEX FOR RUNGE-KUTTA ROUTINE 
COUNTING INDEX 

COUNTING INDEX FOR OUTPUT DETERMINATION 

DUMMY VARIABLE FOR CONVENIENCE 

CONTROL VARIABLE FOR CONVERGENCE ON ANGLE OF 
ATTACK [AOA] 

INCREMENTAL CHANGE IN AOA USING NEWTONS METHOD 
NUMBER OF EQUATIONS OF MOTION 


TIME (SEC) 
FINAL TIME (SEC) 
INCREMENT OF TIME [RUNGE-KUTTA ROUTINE] (SEC) 
TIME PRINT INTERVAL (SEC) 
RADIUS (METERS) 
THETA (RADIANS) 
PHI (RADIANS) 
VELOCITV (M/SEC) 
FLIGHT PATH ANGLE GAMMA (RADIANS) 
PSI (RADIANS) 


DERIVATIVES OF THE ABOVE SIX VARIABLES 
ZEROTH ORDER COEFFICIENT FOR CL EQUATION 
FIRST ORDER COEFFICIENT FOR CL EQUATION 
SECOND ORDER COEFFICIENT FOR CL EQUATION 
ZEROTH ORDER COEFFICIENT FOR CD EQUATION 
FIRST ORDER COEFFICIENT FOR CD EQUATTON 
SECOND ORDER COEFFICIENT FOR CD EQUATION 


DENSITX MODEL EXPONENT (METERS) 
GRAVITY AT EARTH SURFACE (9.806M 2743 
LOCAL GRAVIATIONAL ATTRACTION (М^275) 
REFERENCE ALTITUDE FOR DENSITY MODEL (METERS) 
REFERENCE DENSITY FOR DENSITY MODEL (KG/M^3) 


GRAVITATIONAL CONSTANT 
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NOLO COO O REIN 


ALPHA 
ALFA 
VEC 
EPS 
MASS 
MASSO 
SPI 

R 

V 
RADIUS 
VELOCITY 
GAMMA 
ANGM 
SIGMA 
OINC 
DELTA 
OMEGA 
AS 

AR 

AW 
THR 
AOB 
RHO 

S 

N 

M 
HEAT 
QS 

ED 

CL 
1111. 
DRAG 


ANGLE OF ATTACK 

ANGLE OF ATTACK 

THRUST VECTORING ANGLE 
THRUST VECTORING ANGLE 

MASS OF SPACECRAFT 

INITIAL MASS OF SPACECRAFT 
SPECIFIC IMPULSE OF POWER PLANT 
RADIUS 

VELOCITY 

RADIUS 

VELOCITY 

FLIGHT PATH ANGLE 

ANGULAR MOMENTUM 

ANGLE OF BANK 

INCLINATION 

INTERMEDIATE ANGLE 

LONGITUDE OF THE ASCENDING NODE 
TANGENTIAL ACCELERATION 
RADIAL ACCELERATION 

BINORMAL ACCELERATION 

THRUST 

ANGLE OF BANK 

DENSITY 

REFERENCE AREA OF SPACECRAFT 


(RADIANS) 
(DEG) 
(DEG) 

(RADIANS) 
(KG) 
(KG) 
(SEC) 
(M) 
(M/S) 
(KM) 
(KM/S) 
(DEG) 
(M*2/S) 
(DEG) 
(DEG) 
(DEG) 
(DEG) 
(M/S*2) 
(M/S*2) 
(M/S*2) 
(N) 

(RADIANS) 
(KG/M^3) 
(M^2) 


DENSITY COEFFICIENT FOR HEATING RATE EQUATION 
VELOCITY COEFFICIENT FOR HEATING RATE EQUATION 


STAGNATION POINT HEATING RATE 


(WATTS/M^2) 


DYNAMIC PRESSURE TIMES THE REFERENCE AREA (N*M^2) 


COEFFICIENT OF DRAG 
COEPPRIGIENT OF LIGT 
LIFT ON SPACECRAFT 
DRAG ON SPACECRAFT 


DEFINED FUNCTION FOR ACCEL DUE TO GRAVITY 


G(R)=MU/(R*R) 


(N) 
(N) 


BU "CccccecececceCCCCCCECCECCCCCECCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 


BEEEEEEEEEEEEEEEEECEEEEECEEECCCCECEEECEEECECCCCCCCCCCCCCCCCCCCCCC 


OPENSIEWOT INPUT FILES. .... DATA.DAT — ORBITAL PARMS 
AERO. DAT — VEHICLE PARMS 


OPEN(10,FILE-'data.dat', STATUS=’OLD’ ) 
OPEN(12,FILE-'aero.dat',STATUS-'OLD') 


ci 


ERSCESENGSSESGCCREGECCECCECCCCCCCCCCCCCCCECCCCCCCCCCCCCCCCCCCCCCCCC 


N 


AAAAAAAAAAA 


O0 O €) 


(ህን (OO 


OPEN(13,FILE-'out.dat', STATUS-'NEWT) 


READ(10,1)(X(I),I=1,6),T,TF,H,TPI 
FORMAT(/,/,20Xx,10(/, 20x DEAD 


READ(12,2)(C(I),I=1,6),AOPT,AOB,N,M,S,MASS0,FM,SPI,THR,EPS 
FORMAT(/,/,/,/,16(/,20X,D13.7)) 


NUMBER OF EQUATIONS TO INTERGRATE 


XC) RADIUS 

X(2)....THETA (SPHERICAL COORD. PARAMETER) 
X(3)....PHI (SPHERICAL COORD. PARAMETER) 
X(4)...,VELOCIN 

X(5)....FLIGHT PATH ANGLE GAMMA 

XD... DSi 


INITIALIZATION OF CONSTANTS 


NE=6 

۳ 1-3. 1 + 1 2 0 
MU=3.986012D+14 
G0=9.806D+0 
MASS=MASSO 
VEC=EPS*180./PI 


MAIN PROGRAM 


INDEX=0 
KOUNT=1 


INITIAL CALL TO GET ACCELERATIONS AND CONTROL VARIABLES 
FOR FIRST OUTPUT 


CALL CNTRL(DRAG,LIFT,THR,ALPHA,AOB,X,MASS,S,C,N,M,RHO, EPS) 
CALL ACEL(AS,AR,AW,DRAG,LIFT,THR,EPS,ALPHA,AOB,MASS) 


CALL HDR(X(1),X(4),MASS,THR, VEC, X(5) ) 
CALL WRT(X,T,THR,MASS, ALPHA, AOB, RHO,N,M) 


CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCECCCECCCCCECCECCCECCCCCCCCECER RSs 


tata (Mt dc 


THIS IS THE MAIN BLOCK OF THE PROGRAM CALLING FIRST THE 
CONTROL SUBROUTINE TO DETERMINE THE OUTPUT VARIABLES. NEXT, 
THE ACCELERATION ROUTINE CALCULATES THE RELATIVE 
ACCELERATION TO THE VEHICLE, THEN THE DIFFERENTIAL EQUATIONS 
OF MOTION ARE DEFINED AND A FOURTH ORDER RUNGE-KUTTA ROUTINE 
IS USED TO INTEGRATE THEM. THIS BLOCK IS COMPUTED AT EACH 
TIME INCREMENT H. 


DZ 


eee ee Gee ee er ee eee eeccceccccceccccccccceccccccCcccCCCCCCCCCCCCCCCCC 


100 


200 


e 
e 


QUO Oo @ COO DET) AAAW 


CALL CNTRL(DRAG, LIFT, THR, ALPHA, AOB,X,MASS,S,C,N,M, RHO,EPS) 
CALL ACEL(AS,AR,AW, DRAG, LIFT, THR, EPS, ALPHA, AOB,MASS) 

CALL ORB(X,XDOT,AS,AR,AW,G(X(1))) 

CALL RK4(T,X,XDOT,NE,H, INDEX) 

IF(INDEX .NE. 0) GO TO 100 


DSNERMINE IF FIME TO PRINT TO OUTPUT USING TIME PRINT 
INTERVAL 


IF(KOUNT .LT. IDNINT(TPI/H)) GO TO 300 
CALL WRT(X,T,THR,MASS, ALPHA, AOB,RHO,N,M) 
PRINT STATUS TO TERMINAL 


WRITE(6,7)T,ALPHA, AOB 
FORMAT('4',F8.1,1X,F6.3,1X,F6.3) 


KOUNT=0 
KOUNT=KOUNT+1 


CHANGE MASS OF S/C 
MASS=MASS-THR*H/(SPI*GO) 


TWO TESTS EITHER CAN STOP THE PROGRAM. TEST FOR MASS LESS 
THAN FINAL MASS AND ALSO IF TIME IS GREATER THAN FINAL TIME 


IF(MASS .GT. FM .AND. T .LE. TF) GO TO 100 
PRINT OUTPUT AGAIN IF MASS SWITCH ENDED PROGRAM 
IF(MASS .LE. FM) CALL WRT(X,T,THR,MASS,ALPHA,AOB,RHO,N,M) 


END 


Nr CCCECCEECCOECEUECCECCECCCCCCCCECCCCCCCCCCCCCCCCCCCCCCCCCce 


THIS SUBROUTINE CALCULATES THE CONTROL OF THE FREE VARIABLES 
DURING THE ATMOSHPHERIC FLIGHT. 


ee ee Oe eee eee eee ee ec ec er ececeeececececcecececccececcccccccecceccccc 


SUBROUTINE CNTRL(DRAG, LIFT, THR,ALPHA,AOB,X,MASS,S, 
*C,N,M,RHO,EPS) 

IMPLICIT REALX8(A-H,L-Z) 

DIMENSION X(6),C(6) 


MU=3.986012D+14 
GAMMA=X (5) 


53 


qA OLOO 


NN (000: €Y0 O0 COA ۰ ۳۵ ን 


ር QO ርጋ 


MEN 


R=X (1) 
VEX (4) 


COEFFICIENTS FOR AN EXPONENTIAL DENSITY MODEL 
REFERENCE J MEASE 1976 US STANDARD ATMOSPHERE 
RANGE 50KM TO 120KM 


BETA=1.41D-4 
RHOO-3.0968D-4 
R0=6.435D+6 


G=MU/(R*R) 


RHO=RHOO*EXP(-BETA*(R-RO)) 
QS=.5D+0*RHO*V*V*S 


FOR AEROBANG THE AOB IS SET AND THE ANGLE OF ATTACK IS 
CONTROLLED TO FLY AT A CONSTANT HEATING RATE; THE ALTITUDE 
AND VELOCITY ARE ALLOWED TO FLOAT. 

GUESS ALPHA INITIALLY FOR NEWTON APPROXIMATION ROUTINE 


ALPHA=(C(5)*QS+SQRT(C(5)**2+4.*QS*C(6)*(THR*COS(EPS) 


* -0S*C(4))))/(2.*C(6)*QS) 


USE NEWTON APPROXIMATION METHOD TO CONVERGE ON ALPHA. 
USE A WHILE STRUCTURE FOR THE CONVERGENCE OF ALPHA. 

SET OLDCO AND CHANGE EQUAL TO THE VALUES BELOW TO ENSURE 
THAT THE PROGRAM SWITCHES ARE NOT INITIALLY TRIPPED. 


OLDCO=1.D+10 
CHANGE=1. 


OR, ALLOW ALPHA TO BE A CONSTANT AND OBSERVE THE BEHAVIOR OF 
THE HEATING RATE AND INCLINATION. 


ALPHA = 0.69813 

DO WHILE (ABS(CHANGE) .GT. 1.D-3) 
CL=C(1)+C(2)*ALPHA+C(3)*ALPHA**2 
CD=C(4)+C(5)*ALPHA+C(6)*ALPHA**2 


CALCULATE THE LIFT AND DRAG FORCES FROM THE COEFFICIENTS 
CALCULATED ABOVE. 


LIFT=CL*QS 
DRAG=CD*OS 


THESE THREE EQUATIONS REPRESENT THE FUNCTION AND ITS 
DERIVATIVE WITH RESPECT TO ALPHA FOR NEWTON’S METHOD OF 
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E APPROXIMATION 


fa 

COEFF=(THR*COS(ALPHA+EPS)-DRAG)-MASS*SIN(GAMMA) 

**(G+BETA*(N/M)*V*V) 

COEFP=-THR*SIN(ALPHA+EPS)-(C(5)+2*C(6)*ABS(ALPHA))*QS 
8 

CHANGE=COEFF/COEFP 

ALPHA=ALPHA-CHANGE 

IF (ABS(OLDCO)-ABS(CHANGE).LE.0.) STOP 

OLDCO=CHANGE 

END DO 

RETURN 

END 
С 
C 
NEN UsEECCCCCCCCCCCCCCCCCECCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 
С 
С THIS SUBROUTINE DEFINES THE TANGENTIAL, NORMAL, AND BINORMAL 
C ACCELERATIONS WITH RESPECT TO THE VEHICLE AND THE ORBITAL 
€ PLANE. 
С 
ሀበ... ፡ርሮሮሮርሮሮሮሮሮሮሮሮሮሮሮሮሮሮርርሮርርርርርርርርርርርርርርርርርርርርርርርርርርርርርርርር 
E 


SUBROUTINE ACEL(AS,AR,AW,DRAG,LIFT,THR,EPS,ALPHA,AOB,MASS) 
IMPLICIT REAL*8(A-H,L-Z) 


COMPUTE THE ACCELERATIONS ON THE FLIGHT VEHICLE USING 
THE THREE ACCELERATION EQUATIONS, WHICH INCLUDE THE 
THRUST VECTORING ANGLES. 


OOOO 


AS=(THR*COS(ALPHA+EPS)-DRAG)/MASS 
AW=(LIFT+THR*SIN(ALPHA+EPS))*SIN(AOB)/MASS 
AR=(LIFT+THR*SIN(ALPHA+EPS))*COS(AOB)/MASS 


RETURN 
END 


ፊይ ፤ማ ማ ከ ማዎሉችሽሜከመከከ ከሮሮሮሮማ(.:((ር(ር ር እ ሎ.(..(.ብብ..ብ.......ቂቂከርጀ፡ርርርርርርርርርርርርርርርር 
THIS SUBROUTINE IS THE COLLECTION OF DIFFERENTIAL EQUATIONS 
WHICH DESCRIBE THE MOTION OF THE SPACECRAFT. 


MIMMEBBE7CEEEEEEEĊEECCECECCECECECECCECCCECCECCEECECCCCCOCCCCCCCCCCC 


MOON AA 


SUBROUTINE ORB(X,XDOT,AS,AR,AW,G) 
IMPLICIT REAL*8(A-H,L-Z) 
DIMENSION X(6),XDOT(6) 


XDOT(1)=X(4)*SIN(X(5)) 
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XDOT (2 )=X(4)*COS(X(5) )*COS( x16) C CO ae 
XDOT(3)=X(4)*COS(X(5))*SIN(X(6))/X(1) 
XDOT(4)=AS-G*SIN(X(5)) 
XDOT(5)=(AR-G*COS(X(5))+X(4)*X(4)*C0S(X(5))/X(1))/X(4) 
DUM=TAN(X(3))/X(1) 
XDOT(6)=AW/(COS(X(5))*X(4))-X(4)*COS(X(5))*COS(X(6))*DUM 


RETURN 
END 


CCCCCCCcecccccccccccccceccccecccccececceċceccċccececcċcec c ce c cc ta ደይ ሔዋሎዊሎሮፎ 


E 
C 
C 


THIS SUBROUTINE TAKES CARE OF ALL HARDCOPY OUTPUT 


CCCCCCCccccccccccccccecccccccccccccccccccecccecececcccceclceccc € GG ይይይ 


E 


1 


E 


SUBROUTINE WRT(X,T,THR,MASS,ALPHA, AOB, RHO,N,M) 
IMPLICIT REAL*8(A-H,L-Z) 
DIMENSION X(6),XDOT(6) 


PI=3.14159265359D+0 


RADIUS=X(1)/1000. 

VELOCITY=X(4)/1000. 

GAMMA=X(5)*180./PI 

ALFA=ALPHA*180./PI 

HEAT=RHO**N*X(4)**M*9.652D-5 

SIGMA=AOB*180./PI 

OINC=ACOSD(COS(X(3))*COS(X(6))) 

DELTA=0.0 

IF (OINC .NE. 0.) DELTA=ASIND(TAN(X(3))/TAND(OINC)) 
OMEGA=X(2)*180./PI-DELTA 


WRITE(13,1)T,RADIUS,VELOCITY,MASS,GAMMA,OINC, THR, HEAT, 


*ALFA,SIGMA 


FORMAT(1X,F6.1,1X,F8.3,1X,F6.4,2X,F7.2,1X,F5.3,1X, 


* F6.3,1X,F7.1,1X,E9.3, 1. 57 2, 1 ШЕЛ, 


RETURN 
END 


СССССССССЄССССЄСССССЄСССССССССССССССССССССССССССССССССССССССС LL 


C 
C 
C 
E 


THIS SUBROUTINE ATTACH A HEADER TO THE OUTPUTS TO KEEP TRACK 
OF THE DIFFERENT CHANGES 


ሮርርርርርርርርርርርርርርርርርርርርርርርርሮርርሮርርርርርርርርርርሮሮሮርሮ ርር ርር ርር ርር ር 6 


С 


SUBROUTINE HDR(R,V,MASS, THR, VEC, GAMMA) 
IMPLICIT REAL*8(A-H,L-Z) 
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DO 1 I=13,14 
1 (7.EO.14)C0 TO ı 


WRITE(I,*) 
WRITE(I,*)' SELECTED INITIAL INPUT DATA: ’ 
WRITE(I,*) 
WRITE(I,2)R,V,MASS, THR, VEC, GAMMA 
WRITE(I,*) 
WRITE(I,*) 
1 CONTINUE 
2 FORMAT(’ INITIAL RADIUS (METERS) TIO 2 INITIAL 
*VELOCITY (METERS/SEC) I o INITIAL MASS (KG) 
* S7 INITIAL THRUST (NEWTONS) ' 
*,F9.2,/,' THRUST VECTOR ANGLE (DEG) : 
*,F6.2,/,' INITIAL FLIGHT PATH ANGLE (RADS) ',Е5.2) 
е 
WRITE(13,*)' TIME RADIUS VELOCITY MASS GAMMA INCLI 
* THRUST QDOT ALPHA AOB' 
WRITE(13,*)’ (SEC) (KM) (KM/SEC) (KG) (DEG) (DEG) 
* (N) (0/М2©) (DEG) (DEG)' 
G 
RETURN 
END 
С 


IMMIMIMME EEEEECCEEEECCCECCCEECECECCECECCCCCECCCCCCccccccccccccccccccc 
E 


С THIS SUBROUTINE IS А FOURTH ORDER RUNGE-KUTTA ROUTINE TAKEN 
@ PROM PROF: ROSS’S SUBROUTINE. 
E 


Mune ee GCCCCCCCCCCCCCCECCececccceccccccccccccccccccccccccccccccce 
С 

SUBROUTINE RK4(T,X,XDOT,NE,H, INDEX) 

IMPLICIT REAL*8(A-H,L-Z) 

INTEGER INDEX,I 

DIMENSION X(6),XDOT(6),SAVED(6),SAVEX(6) 


INDEX=INDEX+1 

CORO (1,2,3,4), INDEX 
1 DO 10 I=1,NE 

SAVEX(I)=X(I) 

SAVED(I)=XDOT(I) 


10 X(I)=SAVEX(I)+.5DO*H*XDOT(I) 
T=T+.5D0*H 
RETURN 

E 

2 DO 20 I=1,NE 
SAVED(I)=SAVED(I)+2.DO*XDOT(I) 

20 X(I)=SAVEX(I)+.5DO*H*XDOT(I) 
RETURN 

С 

3 DO 30 I=1,NE 


E 


30 


e eQ 


SAVED(I)=SAVED(I)+2.DO*XDOT(I) 
X(I)=SAVEX(I)+H*XDOT(I) 
T=T+.5D0*H 

RETURN 


DO 40 I=1,NE 
X(I)=SAVEX(I)+H/6.DO*(SAVED(I)+XDOT(I)) 
INDEX=0 

RETURN 

END 
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SAMPLE AERO.DAT INPUT FILE 


This file gives the aerodynamic data used for the MRRV. 


C1-6,AOPT,N,M,S,MASS,FM,SPI,THRUST 


COEFFICENT Cl -1.000000D-02 
COEFFICENT C2 2.860000D-01 
COEFFICENT C3 1.313000D+00 
COEFFICENT C4 4.700000D-02 
COEFFICENT C5 -4.470000D-01 
COEFFICENT C6 2.040000D+00 
AOA FOR CL/CD MAX 1.823870D-01 
AOB FOR BANG MANU 1.570796D+00 
DENSITY EXPONENT 0.500000D+00 
VELOCITY EXPONENT 3.150000D+00 
REFERNCE AREA 1.169800D+01 
INITIAL MASS (KG) 4.898000D+03 
FINAL MASS (KG) 4.800000D+03 
SPECIFIC IMPULSE (S) 2.950000D+02 
THRUST (N) 1.467900D+04 


THRUST ANGLE (RAD) 0.261800D+00 
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SAMPLE DATA.DAT INPUT FILE 


This file gives the orbital and time printing information 
the simulation program. 


INPUT ORBITAL PARAMETERS 


RADIUS METERS 6.4450000D+06 
THETA RADIANS 0.0000000D+00 
PINE RADIANS 0.0000000D+03 
VELOCITY M/S 7.7100000D+03 
FLIGHT PATH ANGLE 0.0000000D+00 
PSI RADIANS 0.0000000D+00 
BEGIN TIME SEC 0.0000000D+00 
END TIME SEC 4.0000000D+03 


TIME INTERVAL SEC 5.0000000D-03 
PRINT TIME INTERVAL 1.0000000D+00 
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SAMPLE OUT.DAT OUTPUT FILE 


SELECTED INITIAL INPUT DATA: 


NITIAL RADIUS (METERS) 6445000.00 
NITIAL VELOCITY (METERS/SEC) 7 7 188۶: 0 
NITIAL MASS (KG) 4898.00 
NITIAL THRUST (NEWTONS ) 14679.00 
HRUST VECTOR ANGLE (DEG) 15.00 
NITIAL FLIGHT PATH ANGLE (RADS) 0.00 
TIME RADIUS VELOCITY MASS GAMMA INCLI THRUST 
(SEC) (KM) (KM/SEC) (KG) (DEG) (DEG) (N) 
0.0 6445.000 7.7100 4898.00 0.000 0.000 14679.0 
1.0 6445.000 7.7100 4892.95 -.003 0.037 14679.0 
2.0 6444.999 7.7099 4887.88 -.006 0.074 14679.0 
3.0 6444.998 7.7097 4882.80 -.008 0.112 14679.0 
4.0 6444.997 7.7095 4877.73 -.011 0.150 14679.0 
5.0 6444.995 7.7092 4872.65 -.014 0.189 14679.0 
6.0 6444.993 7.7088 4867.58 -.017 0.228 14679.0 
7.0 6444.991 7.7084 4862.50 -.019 0.268 14679.0 
8.0 6444.988 7.7079 4857.43 -.022 0.309 14679.0 
9.0 6444.985 7.7074 4852.36 -.025 0.349 14679.0 
10.0 6444.981 7.7068 4847.28 -.028 0.391 14679.0 
11.0 6444.977 7.7061 4842.21 -.031 0.433 14679.0 
12.0 6444.973 7.7053 4837.13 -.034 0.475 14679.0 
13.0 6444.968 7.7045 4832.06 -.036 0.518 14679.0 
14.0 6444.963 7.7037 4826.98 -.039 0.561 14679.0 
15.0 6444.958 7.7027 4821.91 -.042 0.605 14679.0 
16.0 6444.952 7.7017 4816.84 -.045 0.649 14679.0 
17.0 6444.946 7.7006 4811.76 -.048 0.694 14679.0 
18.0 6444.939 7.6995 4806.69 -.051 0.739 14679.0 
19.0 6444.932 7.6983 4801.61 -.054 0.785 14679.0 
19.3 6444.930 7.6979 4799.99 -.055 0.800 14679.0 
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5 5۵ 5۵ 5۵ ር ОО СОО ОСОО О ОСОО ОСОО Оо Оо ዕር ር oo00 © 


ODOT 
(J/M25S) 
.147E+07 
.147E+07 
.147E+07 
21475807 
.147E+07 
.147E+07 
147407 
.147E+07 
.147E+07 
.147E+07 
.147E+07 
.147E+07 
.147E+07 
.147E+07 
.147E+07 
.147E+07 
.147E+07 
.147E+07 
.147E+07 
.147E+07 
.147E+07 


ALPHA 
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